Download NHD layers for Klamath Basin
layers <- st_layers(here::here('data-raw', 'shapefiles', 'NHDPLUS_H_1801_HU4_GDB.gdb'))
nhd_fcode <- st_read(here::here('data-raw', 'shapefiles', 'NHDPLUS_H_1801_HU4_GDB.gdb'), layer = 'NHDFCode') |>
janitor::clean_names()
## Reading layer `NHDFCode' from data source
## `/Users/maddeerubenson/Documents/git/rivermmile/data-raw/shapefiles/NHDPLUS_H_1801_HU4_GDB.gdb'
## using driver `OpenFileGDB'
nhd_flowline_vaa<- st_read(here::here('data-raw', 'shapefiles', 'NHDPLUS_H_1801_HU4_GDB.gdb'), layer = 'NHDPlusFlowlineVAA') |>
janitor::clean_names()
## Reading layer `NHDPlusFlowlineVAA' from data source
## `/Users/maddeerubenson/Documents/git/rivermmile/data-raw/shapefiles/NHDPLUS_H_1801_HU4_GDB.gdb'
## using driver `OpenFileGDB'
nhd_flowline_raw <- st_read(here::here('data-raw', 'shapefiles', 'NHDPLUS_H_1801_HU4_GDB.gdb'), layer = 'NHDFlowline') |>
janitor::clean_names()
## Reading layer `NHDFlowline' from data source
## `/Users/maddeerubenson/Documents/git/rivermmile/data-raw/shapefiles/NHDPLUS_H_1801_HU4_GDB.gdb'
## using driver `OpenFileGDB'
## Simple feature collection with 341762 features and 18 fields
## Geometry type: MULTILINESTRING
## Dimension: XYZM
## Bounding box: xmin: -124.4217 ymin: 38.29819 xmax: -120.5703 ymax: 43.33625
## z_range: zmin: 0 zmax: 1612
## m_range: mmin: 0 mmax: 100
## Geodetic CRS: NAD83
huc_12 <- st_read(here::here('data-raw', 'shapefiles', 'NHDPLUS_H_1801_HU4_GDB.gdb'), layer = 'WBDHU12') |>
janitor::clean_names() |>
st_transform(crs = "+proj=longlat +datum=WGS84") |>
filter(substr(huc12, 1, 6) == "180102")
## Reading layer `WBDHU12' from data source
## `/Users/maddeerubenson/Documents/git/rivermmile/data-raw/shapefiles/NHDPLUS_H_1801_HU4_GDB.gdb'
## using driver `OpenFileGDB'
## Simple feature collection with 705 features and 21 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -124.5351 ymin: 38.25486 xmax: -120.5637 ymax: 43.34273
## Geodetic CRS: NAD83
huc_10 <- st_read(here::here('data-raw', 'shapefiles', 'NHDPLUS_H_1801_HU4_GDB.gdb'), layer = 'WBDHU10') |>
janitor::clean_names() |>
st_transform(crs = "+proj=longlat +datum=WGS84") |>
filter(substr(huc10, 1, 6) == "180102")
## Reading layer `WBDHU10' from data source
## `/Users/maddeerubenson/Documents/git/rivermmile/data-raw/shapefiles/NHDPLUS_H_1801_HU4_GDB.gdb'
## using driver `OpenFileGDB'
## Simple feature collection with 147 features and 16 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -124.5351 ymin: 38.25486 xmax: -120.5637 ymax: 43.34273
## Geodetic CRS: NAD83
huc_8 <- st_read(here::here('data-raw', 'shapefiles', 'NHDPLUS_H_1801_HU4_GDB.gdb'), layer = 'WBDHU8') |>
janitor::clean_names() |>
st_transform(crs = "+proj=longlat +datum=WGS84") |>
filter(substr(huc8, 1, 6) == "180102")
## Reading layer `WBDHU8' from data source
## `/Users/maddeerubenson/Documents/git/rivermmile/data-raw/shapefiles/NHDPLUS_H_1801_HU4_GDB.gdb'
## using driver `OpenFileGDB'
## Simple feature collection with 22 features and 14 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -124.5351 ymin: 38.25486 xmax: -120.5637 ymax: 43.34273
## Geodetic CRS: NAD83
fcode_descr <- c("Stream/River", "Artificial Path",
"Stream/River: Hydrographic Category = Intermittent",
"Stream/River: Hydrographic Category = Perennial", "Stream/River: Hydrographic Category = Ephemeral")
nhd_flowline <- nhd_flowline_raw |>
mutate(rc_huc_6 = substr(reach_code, 1, 6),
rc_huc_8 = substr(reach_code, 1, 8),
rc_huc_10 = substr(reach_code, 1, 10),
rc_huc_12 = substr(reach_code, 1, 12)) |>
filter(rc_huc_6 == "180102") |> # filter to klamath huc6
filter(!is.na(gnis_name)) |> # remove NA streams
left_join(nhd_fcode) |>
filter(description %in% fcode_descr) |>
left_join(nhd_flowline_vaa) |>
glimpse()
## Rows: 31,796
## Columns: 73
## $ permanent_identifier <chr> "137190692", "147512448", "122413487", "8…
## $ f_date <dttm> 2012-02-26 18:16:31, 2012-02-26 18:23:10…
## $ resolution <int> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
## $ gnis_id <chr> "01141109", "01156651", "01657346", "0112…
## $ gnis_name <chr> "Doe Creek", "Long Creek", "Scaath Creek"…
## $ length_km <dbl> 3.35500000, 5.33000000, 0.08269511, 4.451…
## $ reach_code <chr> "18010201000166", "18010202000293", "1801…
## $ flow_dir <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ wb_area_permanent_identifier <chr> NA, NA, "122415943", "84069699", NA, NA, …
## $ f_type <int> 460, 460, 558, 558, 460, 460, 460, 460, 4…
## $ f_code <int> 46003, 46006, 55800, 55800, 46003, 46003,…
## $ main_path <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ in_network <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ visibility_filter <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ shape_length <dbl> 0.0374853923, 0.0606200004, 0.0007596266,…
## $ nhd_plus_id <dbl> 5.00004e+13, 5.00004e+13, 5.00004e+13, 5.…
## $ vpuid <chr> "1801", "1801", "1801", "1801", "1801", "…
## $ enabled <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ rc_huc_6 <chr> "180102", "180102", "180102", "180102", "…
## $ rc_huc_8 <chr> "18010201", "18010202", "18010209", "1801…
## $ rc_huc_10 <chr> "1801020100", "1801020200", "1801020900",…
## $ rc_huc_12 <chr> "180102010001", "180102020002", "18010209…
## $ description <chr> "Stream/River: Hydrographic Category = In…
## $ canal_ditch_type <chr> " ", " ", " ", " ", " ", " ", " ", " ", "…
## $ construction_material <chr> " ", " ", " ", " ", " ", " ", " ", " ", "…
## $ hydrographic_category <chr> "Intermittent", "Perennial", " ", " ", "I…
## $ inundation_control_status <chr> " ", " ", " ", " ", " ", " ", " ", " ", "…
## $ operational_status <chr> " ", " ", " ", " ", " ", " ", " ", " ", "…
## $ pipeline_type <chr> " ", " ", " ", " ", " ", " ", " ", " ", "…
## $ positional_accuracy <chr> NA, NA, " ", " ", NA, NA, NA, NA, NA, NA,…
## $ relationship_to_surface <chr> " ", " ", " ", " ", " ", " ", " ", " ", "…
## $ reservoir_type <chr> " ", " ", " ", " ", " ", " ", " ", " ", "…
## $ stage <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
## $ stream_leve <int> 4, 4, 5, 4, 4, 4, 7, 7, 5, 5, 5, 7, 4, 6,…
## $ stream_orde <int> 2, 4, 2, 6, 2, 2, 3, 3, 2, 3, 2, 3, 3, 1,…
## $ stream_calc <int> 2, 4, 2, 6, 2, 2, 3, 3, 2, 3, 2, 3, 3, 1,…
## $ from_node <dbl> 5.00004e+13, 5.00004e+13, 5.00004e+13, 5.…
## $ to_node <dbl> 5.00004e+13, 5.00004e+13, 5.00004e+13, 5.…
## $ hydro_seq <dbl> 5.00004e+13, 5.00004e+13, 5.00004e+13, 5.…
## $ level_path_i <dbl> 5.00004e+13, 5.00004e+13, 5.00004e+13, 5.…
## $ path_length <dbl> 0.00000, 16.63400, 13.92796, 0.01400, 17.…
## $ terminal_pa <dbl> 5.00004e+13, 5.00004e+13, 5.00004e+13, 5.…
## $ arbolate_su <dbl> 8.332000, 73.236400, 3.444947, 10043.2644…
## $ divergence <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ start_flag <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ terminal_fl <int> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ up_level_pat <dbl> 5.00004e+13, 5.00004e+13, 5.00004e+13, 5.…
## $ up_hydro_seq <dbl> 5.00004e+13, 5.00004e+13, 5.00004e+13, 5.…
## $ dn_level <int> -1, 4, 4, 4, 4, 4, 7, 6, 5, 5, 4, 7, 4, 5…
## $ dn_level_pat <dbl> 0.00000e+00, 5.00004e+13, 5.00004e+13, 5.…
## $ dn_hydro_seq <dbl> 0.00000e+00, 5.00004e+13, 5.00004e+13, 5.…
## $ dn_minor_hyd <dbl> 0.00000e+00, 5.00004e+13, 0.00000e+00, 5.…
## $ dn_drain_cou <int> 0, 2, 1, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ from_meas <dbl> 0.00000, 0.00000, 0.00000, 0.00000, 0.000…
## $ to_meas <dbl> 59.11259, 100.00000, 5.04551, 40.59847, 1…
## $ rtn_div <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ thinner <int> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
## $ vpu_in <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ vpu_out <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ area_sq_km <dbl> 8.67549996, 9.59290002, 0.00760000, 3.879…
## $ tot_da_sq_km <dbl> 12.8034, 87.4999, 1.6544, 10137.9325, 34.…
## $ div_da_sq_km <dbl> 12.8034, 87.4999, 1.6544, 2510.5202, 34.3…
## $ max_elev_raw <dbl> -9998, -9998, -9998, -9998, -9998, -9998,…
## $ min_elev_raw <dbl> 140977, 155705, 655, 124234, 149933, 1489…
## $ max_elev_smo <dbl> 147172, 165024, 800, 124607, 153047, 1644…
## $ min_elev_smo <dbl> 140977, 155705, 655, 124252, 149977, 1489…
## $ slope <dbl> 0.01846498, 0.01748405, 0.01753429, 0.000…
## $ slope_len_km <dbl> 3.35500000, 5.33000000, 0.08269511, 4.451…
## $ elev_fixed <int> 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 1,…
## $ hw_type <int> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
## $ hw_node_sq_km <dbl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
## $ status_flag <chr> "A", "A", "A", "A", "A", "A", "A", "A", "…
## $ Shape <MULTILINESTRING [°]> MULTILINESTRING ZM ((-121…
ggplot() +
geom_sf(data = hucs, aes(color = 'darkred')) +
geom_sf(data = nhd_flowline) +
theme_minimal() +
theme(legend.position = "none")
# filtering to mainstem and tribs of interest
rivers <- c('Klamath River', 'Trinity River', 'Sprague River',
'Williamson River', 'Wood River', 'Lost River', 'Shasta River',
'Scott River', 'Salmon River', 'Blue Creek', 'Bogus Creek', 'Clear Creek',
'Indian Creek', 'Link River', 'Jenny Creek')
nhd_flowline_filtered <- nhd_flowline |>
filter(gnis_name %in% rivers) |>
st_transform(crs = "+proj=longlat +datum=WGS84") |>
st_zm()
leaflet(data = nhd_flowline_filtered) |>
addTiles() |>
addPolylines(popup = ~paste0("name: ", gnis_name, "<br>",
"id: ", gnis_id))
I decided to process the NHD file in QGIS - the project is saved
within the data-raw/shapefiles folder.
Processing Notes:
Merged the trinity river NHD lines at the Whiskeytown Lake Dam
Removed disjointed NDH layers
Connected the Williamson River through the Klamath Marsh
Clipped lines to not start within rivers/reservoirs
Multiple Indian Creeks - kept Indian Creek that is confluence to Trinity River
# Read in processed
processed_nhd <- sf::read_sf(here::here('data-raw', 'shapefiles', 'cleaned_nhd_layer.shp'))
leaflet(processed_nhd) |>
addTiles() |>
addPolylines(popup = ~paste0("name: ", processed_nhd$gnis_nm, "<br>",
"id: ", processed_nhd$gnis_id))
klamath_rivers <- processed_nhd |>
group_by(gnis_nm) |>
summarise(geometry = st_union(geometry)) |>
st_cast("MULTILINESTRING")
# st_transform(crs = 32610)
leaflet(klamath_rivers) |>
addTiles() |>
addPolylines(popup = ~paste0(gnis_nm))
usethis::use_data(klamath_rivers, overwrite = TRUE)
## ✔ Setting active project to "/Users/maddeerubenson/Documents/git/rivermmile".
## ✔ Saving "klamath_rivers" to "data/klamath_rivers.rda".
## ☐ Document your data (see <https://r-pkgs.org/data.html>).
knitr::knit_exit()